CHINESE JOURNAL 
Apr. 2011 Vol. 28 No. 2 





Article ID:1005-3085(2011)02-0272-07 





The Primal-dual Active Set Method for the 
Complementarity Problem* 


LI Qing-guo, YANG Hai-jian 
(College of Mathematics and Econometrics, Hunan University, Changsha 410082) 


Abstract: In this paper, we deal with the convergence properties of a primal-dual active set method for 
the complementarity problem with T-monotone operators. We prove that the primal-dual 
active set method can be interpreted as a specific semismooth Newton method applied to 
this kind of complementarity problems. The established convergence results and numerical 
tests imply that the iteration number of the method is bounded by the number of the 
unknowns. Finally, numerical results show the efficiency of the proposed method. 
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1 Introduction 


The complementarity problem has many important applications in operation research, eco- 
nomic equilibrium models and in the engineering sciences!!:?!_ For this reason, there is a grow- 
ing interest in finding efficient and robust algorithms for solving the complementarity problem. 
This reflects in an increasing number of proposals of solution schemes for the complementarity 
problem in recent years. In these recent developments, an important role has been played by 
the semismooth methods, i.e., by those methods that attempt to solve the complementarity 
problem by first reformulating it as a semismooth system of equations and then by applying a 
generalized Newton method to solve this system. Another efficient way to solve complementar- 
ity problems is given by the primal-dual active set method!3). Its basic iteration consists of two 
steps. In the first phase, based on a certain criterion the domain is decomposed into active and 
inactive parts. In the second phase, a reduced nonlinear system associated with the inactive set 
is solved. We remark that the above methods are based on the linear complementarity problem. 

In this paper, the semismooth Newton method and the primal-dual active set method for a 
kind of nonlinear complementarity problems are described, respectively. It can be thought of as 
an extension for the linear complementarity problem, which was proposed by Hintermiuller et 
alll. The established convergence results and numerical tests imply that the iteration number 
of the method is bounded by the number of the unknowns. 
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2 Preliminaries 


Let N = {1,2,---,n} be an index set. For J, J C N, we denote Azz as the submatrix of 
the matrix A € R”*” consisting of a;; (i € I, j € J), rz as the subvector of r € R” consisting 
of r; (i € I). Let K C R”, f be an operator from K to R”, and for all v € K be expressed by 
v=vt+v7 with vt = max{v,0}, v™ = min{v,0}. Then the notion of T-monotone can be 
defined as follows. 

Definition 2.1 The operator f is called T-monotone over K, if 


(f(v) — f(w), (v-w)t) >0, Vuwek, (1) 


where (-,-) denotes the inner product. Moreover, if for all v, w € K, (f(v) — f(w), (v-w)t) =0 
is equivalent to (v — w)* = 0, then f is called strictly T-monotone over K. 

In this paper, we consider the following nonlinear complementarity problem: find a solution 

u* € R” such that 

ur, f(u)>0, (u—¢)"f(u) =0, (2) 
where f : R” — R” is a continuously differentiable, coercive and strictly T-monotone operator 
and ¢ € R” is a given vector. It is well known that problem (2) has a unique solution!!. 

Lemma 2.1 Let f be a continuous T-monotone operator over K and I, J C N satisfying 
J = N\I. For any vectors y, z € K, if yr = zy and yy > zy, then fr(y) < fr(z). 

Proof Let Î = {i € 1: fily) > filz)} and J = N\f. Without loss of generality, let 
[= {1,2,---,k} and J= {k+1,---,n}. Suppose that Î is not empty. Denote wı = {21 + 
61,°°: Zk + Ôk, Ze41;°°* ; Zn} and w2 = y, where 6; is a positive integer for all i € Î. If 6; is 
small enough, then f;(wi) < fi(we) for alli € I, since the operator f is continuous. Moreover, 


note that (wı — w2)} = 0, we have 


k 
0 < (f(w1) — f(we), (wi — we)t) = Na (fi(wr) — fi(we)) < 0, 
i=l 


which is a contradiction. Hence Î = Ø, which means that fr(y) < fr (z). 
Lemma 2.2 Let f be a continuous strictly T-monotone operator over K and I, J C N 
satisfying J = N\I. For any vectors y, z € K, if yr < zy and fj(y) < f(z), then y < z. 
Proof Let Î = {i € N : yi < zi} and J = N\Î. Without loss of generality, let Î = 
{1,2,---,k} and J= N\i. Suppose that J is not empty, then by the definition of 7, we have 
that 7 c I and J c J. Moreover 


Fiu) < falze) Y> zy yp S 2. (3) 


Hence f5(zj,¥3) < filur ys) = filu) < fil) = Filzi zs), where the first inequality comes 
from (3) and Lemma 2.1. Denote w; = (2;,yj) and w2 = z by the definition of T-monotone 


operator, we have 
0 < (f(w1) — (we), (wi — we)*) = (F5(wi) — F7(we), (yz — z)" <0, 


where the second inequality comes from f ;(w1) < f;(w2) and (y;—z;)t > 0. Hence (fj(wi)— 
f5j(we), (yz — 23)*) = 0, which means y; — z; < 0, since f is a strictly T-monotone operator. 
This is a contradiction to (3), and it means that y < z. 
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3 Semismooth Newton method 


Definition 3.15} Let X and Y be two Banach spaces and L(X,Y) be denoted by the set 
of all linear operators on X into Y. The mapping F : D C X — Y is called slantly differentiable 
in open subset U C D if there exists a family of mappings G: U — L(X,Y ) such that 


lim ne ||F(£ + h) — F(x) — G(x + h)h|| = 0, (4) 
for every x E U. 

We refer to G as a slanting function for F in U. Note that G is not required to be unique 
to be a slanting function for F in U. As shown in [4], the max-function v —> max{0, v} from 
R” — R” is slantly differentiable with a slanting function given by the diagonal matrix Gm(v) 
with diagonal elements 


vü t= 
ms 0, if v; <0. 


In the following, we restrict the Banach space to R”. Then the semismooth Newton method 
is as follows. 

Algorithm 1(semismooth Newton method) 

Step 1 Given an initial vector z? € R”, and set k := 0. 

Step 2 Do the following sequence of vectors 


ght) = gk — gk, (6) 


where ĝz* is the solution of the linear system 
G(z*)z = F(x"), (7) 


with G(x) denoting a slanting function of F at x. 
Step 3 If 2*+! = x*, then stop. Otherwise, set k := k + 1 and return to Step 2. 


In the following, we consider the equivalent form of problem (2) 


f(u) — 


(8) 
ud 9, A> 0, (à,u — ¢) =9, 


Note that the complementarity system given by the second line in (8) can equivalently be 
expressed as 


C(u,A)=0, where C(u,r) = ’— max {0, A+ c(¢—-u)}, (9) 


for each c > 0. Here the max-operation is understood to be componentwise. 


Consequently, (8) is equivalent to 


flu) - 


(10) 
C(u, A) =0 
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Define F : R” x R” — R” x R” by 


F(u,d) = MA , (1) 


A — max{0,\ + c(¢ — u)} 


and note that (10) is equivalent to F(u, A) = 0. As shown in [3], F(u, A) is slantly differentiable 
with the particular choice (5) for the slanting function of the max-function. 

Let (u*,A*) denote the unique solution of (8) and z? = (u°,\°) the initial value of the 
iteration. Similarly to [3], we deduce the following result. For the sake of brevity, the proof is 
omitted. 

Theorem 3.1 Suppose that z* = (u*,A*) is the solution to (10). Then Algorithm 1 


converges superlinearly to x*, provided that ||z° — x*|| is sufficiently small. 


4 The primal-dual active set method 


The primal-dual active set method is based on using (9) as a prediction strategy; i.e., given 
a current primal-dual pair (u, A), the choice for the next inactive and active sets is given by 


T={ieN:r+ce-u)i <0}, A={iEN : M+elo-— u) > 0}. 


This leads to the following primal-dual active set method for problem (2). 
Algorithm 2 (the primal-dual active set method) 
Step 1 Initialize u? and A°. Set k := 0. 
Step 2 Determine the active and inactive sets by 


Ap = {iE N : AF +c(¢—u*); > 0}, (12) 
Th={ieN : X +e- u") <0}. (13) 
Step 3 Let uft! and \**! be the solution of the nonlinear system 
fukt) — Akti = 9, 
u®tl1— gd on Ax, (14) 
Ati = 0 on Ze: 


Step 4 Stop, or set k = k + 1 and return to Step 2. 

We can see that the nonlinear system (14) is the realization of the following nonlinear system 
fa, (u**+?) = 0, 
ut = $a, Agi =0, (15) 
Ne = JA, (utt!) 

Lemma 4.1 Let f be a strictly T-monotone operator. Then (15) has a unique solution. 

Proof Suppose u and v are the solutions of (15). Then we have fr(u) = fr(v), uy = vj = 


J and 
(f(v) — f(u), w- ut) = fre) — frlu), (vr — ur)t) = 0, 
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which means (vr — ur)* = 0, since f is a strictly T-monotone operator. That implies vr < uy. 
Similarly, we have vy > uz. Hence u = v. 

Theorem 4.1 Let f be a continuous strictly T-monotone operator and u* € R” be the 
unique solution of problem (2). Then the sequence {u*} given by Algorithm 2 satisfies 

(i) uf < utt! for all k > 1 and @< u® for all k > 2; 

(ii) u* —> u* at most n iterations. 

Proof Firstly, we show that u* < uft! for all k > 1. Observe that for all k > 1 the 
complementarity property 


AE =0 or uk = Qi, VieN, k<l (16) 


holds. For i € Ap, we have AF + a u*), > 0, and hence by (16) either AF = 0, which implies 
uk < $j, or A£ > 0, which implies uk = ġi. In particular, for i € A, we have 


M>O or uk < 4; (17) 
at each iteration k > 1. In a similar fashion, for i € Z, we obtain 
M<O or uf > di. (18) 


Note that, for i € Tp, fi(u**1) — fi(u®) = Aft? — AF = 0 - AF > 0, where the first equality 
comes from the first line in (12), and the inequality comes from (18). And uk< Qi < uft! 
holds for i € Ap. Hence, we have u¥ < uft! for all k > 1 by Lemma 2.2. 

Next we show that u* is feasible for all k > 2. Due to the monotonicity of uf it suffices to 
show that u? > ġ. Let V = {i € N : ul < ġ}. For i € V we have Xj = 0 by (16), and hence 
Al + c(ġ — ut); > 0, which implies i € Ai. Since u? = ġ on A, and u? > u’, it follows that 

u? > $. 

Turning to the feasibility e A*, assume that for a pair of (k,i), k > 1, we have AX < 0. 
Then it is necessary i € Ag1, UF = Q; and Ak +c(ġ— u"); < 0. It follows that i € Tg, Neth = =0 
and akti +c(d— u"); < 0, since ae > Qi, k > 1. Consequently, i € Zg,, and by induction 
i € Tp for all k >k+1. This means that Zp C Zx41. 

We now prove that the method terminates in at most n iterations. Since T, can not contain 
more than n elements, it follows from Zy C Tp+ı that there is an index kọ € N such that 

ko = Lko+1- Using (18) with k = ko + 1, we obtain 


ko+1 _ „ko+1 
UTi T Uz, et > Pko+1- 


Similarly, it follows from (17) that 


ath = sO eee “2 
Since (uL +! — OVA = 0 and Az,, = 0 in view of (14), respectively, we see that u*t? is a 
solution of problem (2), which is unique from Lemma 4.1. 

Remark 4.1 We can choose u? = ¢ and à? = f(¢) as initialize iterate vectors. From the 
proof of Theorem 4.1, Zk = Tk+1 can be seen as a terminal condition of Algorithm 2. We can also 
conclude that Algorithm 1 based on (7) is equivalent to Algorithm 2 for the complementarity 
problem with T-monotone operator. For the sake of brevity, the proof is omitted. 
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5 Numerical experiments 


Let 2 = (0, 1) x (0, 1) and consider the following nonlinear complementarity problem!l: find 
u € K, such that 


u>0, —Autf(u,z,y)>0, u?(-— Aut f(u,z,y)) =, (19) 


where K = {u : u > 0}, f(u, x,y) = u/(1 + u) +102 +y — 8. We discretize (19) by using the 
standard five-point finite difference method on a uniform grid: h = 1/(n + 1), where n denotes 
the number of mesh nodes in x— or y—directions (N = n? is the total number of unknowns). 
Then (19) reduces to the complementarity problem with the T-monotone operator!4), 

We have conducted the following experiments: comparing the Algorithm 2 with the project 
successive overrelaxation (denoted by PSOR) method!®! and the multiplicative Schwarz (de- 
noted by Schwarz) method!l. In the first test, we consider the numerical solution of (19). To 
be precise, we use PSOR to solve (19) with N = 400. Numerical results are listed in Table 1. 
From the table, we see that the optimal factor is about w = 1.8. 


Table 1: Affect of relaxation parameter w in PSOR method 




















w iter. cpu w iter. cpu 
0.1 7111 13.015 0.3 2396 4.688 
0.5 1345 2.859 0.7 868 1.984 
0.9 591 1.515 1.1 408 1.171 
1.3 275 0.921 1.5 172 0.703 
1.7 79 0.578 1.75 54 0.515 
1.8 52 0.500 1.85 57 0.500 
1.9 81 0.547 1.99 670 1.656 








We compare the three algorithms from the view of iteration numbers and execution times. 
In Schwarz, we partition N = Nı U N3 into two equal parts with the overlapping size O(+), 
and the corresponding subproblems are solved by PSOR with the relaxation parameter w = 1.8. 
In Algorithm 2, the parameter c = 1 and the nonlinear systems are solved by the Gauss-Seidel 
method. We mainly consider the affect of dimension (denoted by N) to the performance of 
each algorithm. In each experiment, K = {u € RY : u > 0} and u? = 0. The tolerance in the 
subproblems of Schwarz was chosen to be equal to 1074 in the ||- |jz-norm, while in the outer 
iterative processes of Schwarz was chosen to be equal to 10~® in the ||- ||2-norm. The tolerance 
in the nonlinear system of Algorithm 2 was chosen to be equal to 10~® in the ||- |]2-norm. 

Table 2 gives the history of iterations and execution times for the above iterative methods. 
From the view of iterations, we can easily see that Algorithm 2 is better than PSOR and Schwarz 
in iteration numbers. Moreover, the behavior of Schwarz is better than that of PSOR. From 
the view of execution times, the domain decomposition approach does not seem reasonable if 
there is no multiprocessor system to be used. And the behavior of Algorithm 2 is much better 
than any other solution methods regardless of the dimension being large or small. 
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Table 2: Comparison of iteration numbers and execution times (iter./cpu) 














N PSOR Schwarz Algorithm 2 
100 34/0 18/0.046 4/0.015 
400 46/0.464 33/1.471 5/0.593 
900 103/4.625 52/11.359 7/5.593 
1600 221 /26.687 74/53.609 10/31.859 

2500 368/102.140 103/225.781 12/120.125 
3600 548/325.234 134/752.250 14/450.609 


From the above numerical results, we can conclude that the proposed methods are effective 


for the complementarity problem with the T-monotone operator. 
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ARB Fp (5) RB AY aT Bee EA Se BSE 
ERE, aie 
GH KARE Si Ae BE, Ke 410082) 


BE: EACH BAN AE EN TREA TSE T EA BT RIE ER. RAX 
ARREA OR APES HAN, YAS AY DAE RR RE GE AE. AC TAL BL 
BA T PES RR BOS EE AR BT. ES, TER AS AT HE. 
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